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Abstract 

The dynamics of solutions of the Einstein-Vlasov system with Bianchi I symmetry is 
i/-) . discussed in the case of massive or massless particles. It is shown that in the case of massive 
particles the solutions are asymptotic to isotropic dust solutions at late times. The initial 
singularity is more difficult to analyse. It is shown that the asymptotic behaviour there 
■ must be one of a small set of possibilities but it is not clear whether all of these possibilities 
are realized. One solution is exhibited in the case of massless particles which behaves quite 
differently near the singularity from any Bianchi I solution with perfect fluid as matter 
model. In particular the matter is not dynamically negligeable near the singularity for this 
' solution. 

On 

"o 
cr 

1. Introduction 

The simplest of all anisotropic cosmological models are those of Bianchi type I. They 
are the spacetimes which admit a three-dimensional abelian symmetry group whose orbits 
are spacelike. (For general information on Bianchi models see [1].) Just how simple their 
dynamics is depends significantly on the nature of the matter content of the spacetime. 
For a perfect fluid with a linear equation of state it has been known for a long time how 
to analyse the dynamics [2,3]. For a non-interacting mixture of two fluids with linear 
equations of state the time evolution is also well understood and is asymptotic near the 
singularity and at large times to that of a single fluid [4] . The case of a fluid with nonlinear 
equation of state is discussed in an appendix to the present paper. The dynamics does 
not differ much from the picture in the linear case. When a magnetic field is added to the 
fluid things are already more complicated. In fact, as was shown by Collins [5], a Bianchi 
type I model with fluid and magnetic field resembles a model of the more complicated 
Bianchi type II with fluid alone. It is also interesting to note that models of type VIo with 
perfect fluid and a magnetic field have a dynamical behaviour resembling the notoriously 
complicated 'Mixmaster' behaviour of Bianchi type IX models [6]. Thus changing the 
matter model can have effects on the complexity of the dynamics comparable with those 
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encountered when passing to more general symmetry types. 

A matter model for which the details of the global dynamics of Bianchi type I space- 
times has not previously been studied mathematically is the collisionless gas, described by 
the Vlasov equation. The only general facts which are known are that, with an appropriate 
choice of time orientation, 

(i) the spacetime is future geodesically complete (when maximally extended towards the 
future) 

(ii) there is a crushing singularity in the past where, except in the vacuum case, the 
curvature invariant R a pR a)3 tends to infinity 

These fundamental facts were proved in [7], where it was shown that they hold for any 
Bianchi type other than IX and for a general class of matter models. The aim of this paper 
is to refine (i) and (ii) in the case of Bianchi type I symmetry and matter described by 
the Vlasov equation so as to get more detailed information about the asymptotics of the 
expanding phase and the nature of the initial singularity. An aspect of the situation which 
makes this more difficult than in the case of many other matter models is that for general 
initial data it is not possible to derive an explicit closed system of ordinary differential 
equations which describes the dynamics. This is because certain integrals which occur can 
not be evaluated. In one special case, where massless particles are considered and the initial 
phase space density has the form of the characteristic function of a ball, these integrals 
have been computed by Lukash and Starobinski[8]. However, the explicit expressions they 
obtain are sufficiently complicated that they do not seem to make a rigorous analysis of 
the global dynamics any easier. On the other hand they would probably be useful for 
numerical calculations since they would allow costly numerical evaluation of integrals to 
be avoided. 

The dynamics at late times of the models with massive particles can be described 
precisely. All solutions become isotropic and can be approximated by dust solutions in 
this limit (Theorem 5.4). On the other hand, the results of this paper do not give a 
complete picture of the dynamics near the initial singularity of the spacetimes being stud- 
ied. They merely reduce the possible types of asymptotic behaviour to a small number 
of alternatives. Improving on this is likely to require new techniques. These results leave 
open the possibility that Bianchi I spacetimes with a matter content described by kinetic 
theory may show complicated oscillatory behaviour and thus be very different from those 
with other types of matter content studied up to now. The mechanism which allows for 
this complexity is simply the presence of anisotropy in the pressure which may respond to 
changes in the geometry. It may be that the only reason that the dynamics is so simple 
in the case of a perfect fluid is that this mechanism is excluded by a special symmetry 
assumption (the isotropy of the pressure). The one conclusion which emerges and which 
applies to all solutions considered here is that the ratio of the mean pressure to the energy 
density tends to one third as the singularity is approached. This means that in a certain 
weak sense the dynamics for particles of unit mass is approximated near the singularity 
by that for massless particles. For this reason both cases are often considered together in 
the following, although the main emphasis is on the case m = 1. 

The results will now be summarized. There are broadly speaking two possible types of 
asymptotic behaviour of solutions of the Einstein- Vlasov system with Bianchi I symmetry 
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near the singularity. They will be referred to as convergent and oscillatory. Let Aj denote 
the eigenvalues of the second fundamental form of the homogeneous hyper surfaces. Then 
the mean curvature of the homogeneous hypersurfaces is given by H = Ai + A2 + A3. Let 
Pi = \/ H. In the convergent type the pi tend to limits as the singularity is approached. 
There are three different cases, depending on these limiting values. The first case is 
that where the limiting values are (|, |, |). The well known homogeneous and isotropic 
solutions of the Einstein- Vlasov system [9] are of this type. The second is that the limiting 
values are (0, |, |) or some permutation thereof. The existence of solutions of this kind in 
the case of massless particles is shown in Section 6. These limiting values of the Kasner 
exponents are not realized by any Bianchi type I spacetime when the matter model is 
a perfect fluid (see appendix). The third is that the limiting values satisfy the Kasner 
relation p\ J rp\ J rp\ = ^-- Any solution for which one of the eigenvalues becomes negative 
at some time has this asymptotic behaviour and so there are plenty of examples. This is 
proved in Theorem 5.1. Note that the special case of this result when two of the pi are 
equal is closely related to the homogeneous special case of a result of Rein[10] for plane 
symmetric spacetimes. In the oscillatory type the pi undergo infinitely many oscillations, 
in a sense which will now be specified. There are two cases to be considered, according to 
whether two of the eigenvalues are always equal or not. Consider first the case where two 
eigenvalues are equal and suppose without loss of generality that A2 = A3. Associate to any 
solution a string of symbols (which may be finite or infinite, depending on the solution) 
as follows. Moving backwards from some fixed time add an x to the string each time that 
Ai — A2 changes from being < to being > and add a y each time it changes from being 
> to being < 0. That this makes sense follows from the fact, proved in Section 2, that 
the set of times where Ai = A2 can have no limit point unless this equality holds at all 
times. Thus a finite time interval can only contribute a finite number of symbols. The 
solution is said to undergo infinitely many oscillations if the resulting string of symbols is 
infinite. Similarly, if there is some time where all eigenvalues are distinct then a string of 
symbols is associated to the solution by adding x, y or z each time Ai, A2 or A3 respectively 
becomes strictly larger than the other two eigenvalues. 

Unfortunately it could not be shown whether any oscillatory solutions exist. If they 
did then the behaviour of Bianchi I models with kinetic theory as matter model would be 
much more complicated than in the case of a perfect fluid. If it could be shown that they 
existed the question would remain whether the sequences of symbols they produce have 
some regularity or whether they are chaotic. 

To each type of solution discussed above corresponds a characteristic behaviour of the 
pressures. The solutions considered in the following all have diagonal energy-momentum 
tensors and so three pressures Pi are defined by three diagonal components. The quantities 
Ri = Pi/p must have a sum which converges to unity at the singularity. When the limiting 
values of the pi are or (0, |, |) then the limiting values of the Ri are (|, |, ^) or 

(0, 2, \) respectively. When the sum of the squares of the limiting values of the pi is equal 
to unity then the Ri tend to (0, |, |) or a permutation thereof, unless one pi has limiting 
value zero. In the latter case the Ri tend to (1, 0, 0) or a permutation thereof. 

The paper is organized as follows. In Section 2 some basic facts about the solutions 
are collected. Section 3 is concerned with a simplified system which in some cases models 
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the asymptotic behaviour of the solutions of the original system. Some estimates for the 
pressures are derived in Section 4. Section 5 contains the main results. Section 6 con- 
tains proofs of the existence or non-existence of solutions with certain kinds of asymptotic 
behaviour. 

2. Basic facts 

The Einstein-Vlasov system for particles of mass m > can be written in the following 
form in the case of Bianchi type I symmetry: 



- kijk ij + H 2 = 16np 


(2.1) 


T 0l = 


(2.2) 




(2.3) 


dthj = Hkij - 2kuk\ - &itT %j - Anpg i3 + 47rtrT#y 


(2.4) 


df/dt + 2kydf/dv i = 


(2.5) 


p = J f{t,v k ){m 2 + g rs v r V s ) 1 / 2 (det g) 1 / 2 dv 1 dv 2 dv 3 


(2.6) 


T ij = J f(t,v k )viVj(m 2 + g rs v r v s )~ 1/2 (det g) 1/2 dv 1 dv 2 dv 3 


(2.7) 



Here gij is the induced metric of the homogeneous hypersurfaces, kij is the second fun- 
damental form, / is the phase space density of particles, p is the energy density, T 0i and 
Tij are components of the energy-momentum tensor and H is the mean curvature g^kij. 
With the exception of / all these quantities depend only on the time coordinate t. This 
time coordinate is Gaussian i.e. it is constant on each homogeneous hypersurface and 
defines a parametrization by proper time when restricted to any geodesic normal to these 
hypersurfaces. The spacetime metric is of the form 

ds 2 = -dt 2 + g lJ (t)dx i dx J (2.8) 

In the following it is always assumed, when talking about solutions of (2.1)-(2.7), that the 
function f(t,v) is non-negative and has compact support for each fixed t. It is assumed 
that / is C 1 except in Section 3, where / may be a distribution. The case of primary 
interest here is the case m = 1. Since, however, solutions of (2.1)-(2.7) with m = 1 
resemble solutions with m = close to the singularity, it is useful to also allow the case 
m = from the beginning. 

For a given Bianchi I geometry the Vlasov equation can be solved explicitly. (This 
is not possible for the other Bianchi types. The reason is explained in [11].) The result 
is that if / is expressed in terms of the covariant components Vi then it is independent of 
time. This means that if t is some fixed time and fo(vi) = f(to,Vi) then (2.6) can be 
rewritten as 

P = J fo(vi)(m 2 + g rs v r v s ) 1/2 (detg)- 1 / 2 dv 1 dv 2 dv 3 (2.9) 

and (2.7) can be rewritten in a similar way. The explicit solution allows certain special 
subclasses of solutions of (2.1)-(2.7) to be identified. The first of these will be referred to 
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as reflection-symmetric and is defined by the conditions that 

fo{vi,v 2 ,v 3 ) = fo{-vi,-v 2 ,v 3 ) = fo{v 1 ,-V2,-v 3 ). (2.10) 

and that the initial values of gij and kij are diagonal. Equation (2.10) implies that is 
diagonal. It then follows from (2.3) and (2.4) that gij and kij are always diagonal. The 
second case, which will be referred to as LRS (locally rotationally symmetric) is obtained 
by supplementing (2.10) by the conditions that 

fo(v u V2, v 3 ) = F(( Vl ) 2 + (v 2 ) 2 , v 3 ) (2.11) 

for some function F and that gn = (722, &11 = ^22 initially. Equation (2.11) implies that 
T11 = T22 and it follows from (2.3) and (2.4) that gn = (722 and k\\ = /C22 everywhere. A 
solution will also be called LRS if it satisfies the definition obtained from that just given 
by a permutation of the indices 1, 2, 3. A solution will be called isotropic if 

/o0>i, i>2, vs) = F(( Vl ) 2 + (v 2 f + (vs) 2 ) (2.12) 

and if g^ and k^ are proportional to Sij on the initial hypersurface and hence everywhere. 
If a solution satisfies the conditions on g^ and k^ in the definition of reflection-symmetric, 
LRS or isotropic for all t in some interval, but no assumption is made on / then the solution 
will be said to have reflection-symmetric, LRS or isotropic geometry respectively on that 
interval. It follows from (2.4) that Ty has a corresponding symmetry property. 

In this paper only reflection-symmetric solutions of (2.1)-(2.7) are considered. The 
alternative notation a 2 = gii/gn(to), b 2 = (722/(722(^0), c 2 = (733/(733 (to) is used when it is 
convenient. Yet another form of (2.6) can then be obtained by doing a change of variables 
in (2.9). 

p = J fo(awi,bw2, cw 3 )(m 2 + 5 rs w r w s ) 1 / 2 dwidw2dw3 (2-13) 

The geometric interpretation of the Wi is that they are the components of the momentum 
in an orthonormal frame. Similarly: 

T\ = J f (aw 1 ,bw2,cw 3 )w 2 (m 2 + 5 rs w r w s )~ 1/2 dw 1 dw 2 dw 3 (2.14) 

In order to study the dynamics of the solutions of the system (2.1)-(2.7) in detail it is 
useful to introduce certain dimensionless variables which remain finite at the singularity. 
It follows from (2.1) that H never vanishes except in the case of flat spacetime, which is 
excluded from consideration in the following. By replacing t by —t if necessary, it can be 
arranged that H < everywhere. It will be assumed that this has been done. Define 



k tJ = kij/H (2.15) 

p = p/H 2 (2.16) 

T ld =T l3 lR 2 (2.17) 

r(£) = - / H(t')dt' (2.18) 

Jt 
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In terms of these variables equations (2.1) and (2.4) become: 



- kfkf + 1 = 167TP 



\5)) + 8nfj + Aix{k) - <yj)trf 



(2.19) 
(2.20) 



The following lemma provides some information about the range of r. 

Lemma 2.1 Let a solution of (2.1)-(2.7) with m = or m = 1 be given which is the 
maximal globally hyperbolic development of initial data on the hypersurface t = to- Then 
H{t) is a monotonic function defined on an interval (£i,oo). By translating t it can be 
assumed that t\ = 0. Then lim^o H = — oo and lim^oo H = 0. Moreover, — 3/t < 
H{t) < -1/t. 

Proof That H(t) is monotonic and defined on an interval of the form (£i,oo) with 
hm t ^ tl H = — oo was shown for the case m = 1 in [12]. Essentially the same argu- 
ment applies for m = 0. In the latter case the coefficients of the characteristic system 
are only Lipschitz instead of C 1 but this causes no problems. Now it follows from (2.4) 
that \H 2 < dtH < H 2 . Comparing the solution with the ordinary differential equations 
corresponding to these inequalities then gives the desired estimates (cf. [7]). 

This result implies that the integral defining r diverges as t — > and as t — > oo. Hence the 
solution of (2.19)-(2.20) exists globally in r. 

The solution of (2.19)-(2.20) of course contains only a small part of the information 
of that contained in the solution of (2.1)-(2.7). The former is only a projection of the 
latter. Nevertheless it will be seen that a lot of information about the solution of the full 
equations can be obtained by studying this projection. Consider now the set K of triples 
of real numbers kl/k^/ki which satisfy ^(&1) 2 < 1 and = 1- This is a compact 

subset of R 3 . In fact it is a disc in a plane. A solution of (2.19)-(2.20) defines a point 
of K at each time r. It is on the boundary of K in the plane if / is identically zero and 
in the interior otherwise. This point depending on r, considered as a mapping from R to 
K, will be referred to as the projection of the given solution. The projection of a vacuum 
solution is constant. Let C denote the point (|, |, |). If a solution has isotropic geometry 
on a time interval then its projection lies at the point C during this time. Conversely, 
if its projection lies at C on a given time interval then it can be made to have isotropic 
geometry by a time-independent rescaling of the spatial coordinates. Let L\, Li and L3 
be the subsets of K defined by k 2 = k%, k\ = k\ and k\ = k\ respectively. A solution has 
LRS geometry on a time interval (up to a constant rescaling of the spatial coordinates as 
above) if and only if its projection lies on one of the lines Li, L 2 or L 3 during this time. 
Let Lf denote the open half of Lj which ends at the point with coordinates (—5, §, §) or a 
permutation thereof and let L~ denote the opposite half-line, which ends at the point with 
coordinates (1,0,0) or a permutation thereof. Let and V~ denote these endpoints. 
Let Ai be the open region bounded by L^", and the boundary of K and let A2 and As 
be defined by cyclically permuting (1,2,3) in this definition. Let Bi be the subset of K 
where k\ < 0. 
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The components of the metric satisfy the evolution equations 

dgu/dr = 2k\g u , d/ dr(g u / g J3 ) = 2{k\ - k J J )(g tt / g 33 ) (2.21) 

which imply that gu or their ratios increase or decrease exponentially if certain sign con- 
ditions are satisfied by the k\. There are of course corresponding statements for the scale 
factors a, b and c. Given an initial datum fo for / the quantities p and T\ are determined 
uniquely by the g„ by means of equations (2.13) and (2.14). Thus (2.19), (2.20) and (2.21), 
together with the equation 

d T H = -H(l - Yl-Kp + 47rtrf ) (2.22) 

derived from (2.4) form a closed system of ordinary differential equations, which formally 
determine the quantities k\, gu and H as functions of r in terms of initial data. If the coef- 
ficients of this system were locally Lipschitz, it would follow from the standard uniqueness 
theorem for ordinary differential equations that they determine them uniquely. It will now 
be shown that in fact for m = 1 this dependence is analytic. To do this it is convenient 
to use the expressions for p and T\ of the type (2.9). Analyticity is a consequence of the 
following lemma: 

Lemma 2.2 Let W be a mapping of U x R 3 to R, where U is an open subset of R 3 . 
Suppose that W extends to a C 1 mapping W of U x R 3 to C, where U is an open 
neighbourhood of U in C 3 , and that W(-,y) satisfies the Cauchy-Riemann equations for 
each fixed y <E R s . Finally suppose that each z G W has an open neighbourhood V such 
that the supports of the functions W(z, •) are contained in a common compact subset K 
of R 3 . Then the function F(x) = J R3 W(x,y)dy is analytic. 

Proof It suffices to show that the function F(z) = J R3 W(z, y)dy is complex analytic and 
this is true if F is C 1 and satisfies the Cauchy-Riemann equations [13]. The assumptions 
on the smoothness and support of W justify differentiation under the integral and so the 
Cauchy-Riemann equations for F follow from the Cauchy-Riemann equations for W. 

A consequence of the analyticity of the coefficients in the system of ordinary differential 
equations is that the solutions are analytic. It follows that if a solution with m = 1 has 
LRS or isotropic geometry on some non-empty open time interval then it must have LRS or 
isotropic geometry respectively for all values of r. The same conclusion holds if there is a 
sequence of times having a limit point where the geometry is LRS or isotropic, respectively. 

3. The asymptotic system 

In this section a certain system of ordinary differential equations is introduced and 
the qualitative behaviour of its solutions analysed. This sytem is used later to study the 
asymptotic behaviour of solutions of the system (2.1)-(2.7). This system can be obtained 
formally from (2.1)-(2.7) by replacing the C 1 function f(t,vi,V2,vs) by a measure of the 
form /(£, Vi)5(v2—V2(t))5(v3—V3(t) where 5 is a Dirac measure and taking m = 0. Solutions 
of this system can be interpreted as certain distributional solutions of the Einstein- Vlasov 
system for massless particles. These are intermediate between smooth solutions and the 
even more singular solutions which are in one to one correspondence with dust solutions. 
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(For the correspondence between dust and distributional solutions of the Vlasov equation 
see [14].) The mathematical results which will now be derived are independent of this 
interpretation. 

The system of ODE's to be considered is the special case of the equations (2.19) and 
(2.20) obtained by setting Tj 1 = T| = and Tf = p. Note that, in contrast to the general 
case of (2.19) and (2.20), these specialized equations suffice to determine all unknowns 
occurring in them from initial data. The explicit form of (2.20) in this case is: 



d T k\ = 


Sirpkl 


(3.1) 


df k<2 ~~ 


— 87T/3/c| 


(3.2) 


d T kl = 


-SivpCkl - 1) 


(3.3) 



If initial data are chosen at some time which satisfy the condition ^i(^l) = 1 then the 
solution also satisfies it. Only solutions with this property are considered here. It follows 
from (3.1) and (3.2) that d T {k\lk\) = whenever k\ ^ 0. Let r be the constant value of 
kl/k\. Then 

k\ = rk\ (3.4) 
fcf = 1 - (1 + r )kl (3.5) 

Substituting (2.19), (3.4) and (3.5) into (3.1) gives: 

d T %l = {kl) 2 [(l + r + r 2 )k\ - (1 + r)\ (3.6) 

Proposition 3.1 Let (kl/k^/kl) be a solution of (3.1)-(3.3) satisfying = 1 an d 

T, t (H) 2 < !• Define r = kl/k\ whenever h\ ^ 0. Then: 

(i) if k\ is zero initially it is always zero 

(ii) if kl is initially (and hence always) non-zero then r is constant 

(hi) when kl ^ it is a monotonic function with linv^-oo kl = (1 + r)/(l + r + r 2 ) and 
linv^oo kl = 

(iv) in that case lim T ^_ 00 (/c3 — k%) = —r(r + 2)/(l + r + r 2 ) 

Proof Statement (i) follows from (3.1). Statement (ii) has been demonstrated above. 
Statement (iii) follows from (3.6). The last conclusion is then an immediate consequence 
of the definitions. 

Of course analogous statements hold if k\ and k\ are interchanged since these two quantities 
occur symmetrically in the hypotheses. 

4. Pressure estimates 

The results of this section are all variations on the theme that if the spacetime is 
expanding in a certain direction then the pressure in that direction tends to decrease. It 
is assumed throughout that the geometry is reflection-symmetric. A solution of (2.1)-(2.7) 
with m = satisfies p = trT. This condition never holds when m = 1 but the next lemma 
shows that it does hold asymptotically. 
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Lemma 4.1 Suppose that some solution of (2.1)-(2.7) with m = 1 is defined on an interval 
(—00, Ti), with / not identically zero. Then lim,-^-^ trT/p = 1. 

Proof It was shown in [7] that lim T ^_ 00 p = 00. Thus the result will follow if it can be 
shown that p — > 00 implies trT/p — > 1. To do this choose some radius L > and write the 
p = P1+P2 and trT = (trT) 1 + (trT) 2 where the first summand is the integral over the region 
\w\ < L of the integrand in (2.13) or (2.14) respectively and the second is the integral over 
the complementary region. Using the fact that(l + a; 2 ) 1 / 2 — x 2 / (1 + x 2 ) 1 / 2 = (1 + x 2 ) -1 / 2 
it can be seen that p 2 — (trT) 2 < (1 + L 2 )~ 1 p 2 . Hence 

(tiT)a > Y^P2 (4.1) 

and 

trT > trT 2 > (p _ (4tt/3)L 3 (1 + L 2 ) 1/2 |l/o||oo) (4.2) 

By choosing L sufficiently large the quantity L 2 /(l + L 2 ) can be made as close to unity as 
desired. For fixed L the quantity in brackets on the right hand side of (4.2) approaches p 
as p becomes large. This suffices to give the conclusion of the lemma. 

For a given initial datum f the equation (2.14) defines the pressures T\ as functions 
of a, b and c. The following results concern the qualitative behaviour of these functions. 

Lemma 4.2 If f is not identically zero and a < C"min{l, 6, c} for some constant C > 
there exists a constant C > such that Tj 1 > Ca -2 6- 1 c -1 and T 2 jT\ < C(a/b) 4 / 3 . In the 
case m = the conclusion holds under the weaker hypothesis that a < C min{6, c} 
Proof Let p be a point of R 3 where foj^O whose first coordinate w± is non-zero. Let 5 
be a positive number such that / is bounded below by some positive constant r] on the 
closed cube W of side 26 centred at p and such that wi does not vanish anywhere on 
this cube. Consider now the image W of the cube W under the mapping (wi, W2, W3) 1— > 
(a~ 1 wi, b~ 1 W2, c~ 1 ws). On W the functions W2/W1 and ws/wi are bounded. Under the 
assumptions of the lemma they are bounded by the same constant on W' . It follows that 
|^i |/(1 + l^l 2 ) 1 / 2 is bounded below on any such cube by a positive constant which is 
independent of a, b and c which satisfy the hypotheses of the lemma. The integral defining 
Ti can be bounded from below by the integral of the same quantity over W. It follows 
that 

Tl > Crja-^c- 1 (4.3) 

and this proves the first part of the lemma. To get a lower bound for T\IT^ the domain 
of integration in the definition of these two quantities will be divided into the regions 
\w2\ > R\wi\ and \w2\ < R\w±\, where R is a positive number which will be specified later. 
Corresponding to this decomposition of the domain of integration, there are decompositions 
Tl = + T\ 2 and T| = T 2 \ + T| 2 . The volume of the region where \w 2 \ > #K| and 
f(awi,bw2,cws) 7^ can be bounded by an expression of the form CR~ 1 c~ 1 b~ 2 and so 
T 2 2 ! < CR- X c- X b- Z . On the other hand, T| 2 < R 2 T\ 2 . Thus 

T 2 < CR-^b' 3 + R 2 Tl < (CR-^a/b) 2 + R 2 )Tl (4.4) 
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where in the last step (4.3) has been used. Choosing R = (a/6) 2/3 gives T| < C(a/6) 4/3 T 1 1 
and this proves the result for Tf/T 1 . 

Lemma 4.3 Suppose that some solution of (2.1)-(2.7) is defined on the interval (— oo, ti), 
with / not identically zero. If k\ — k\ > A and k\ — k 3 > A on this interval for some A > 
then T//T 1 < Ce 4Ar / 3 for % = 2, 3 

Proof It suffices to note that under the assumptions of this lemma there will be a time 
interval (— oo, t 2 ) where the hypotheses of Lemma 4.2 hold, so that (4.4) can be applied. 

The time derivatives of the quantities T\ cannot in general be expressed in terms of 
the dimensionless quantities (2.9)-(2.11) so as to get a closed system of ordinary differential 
equations. However they can be estimated in terms of these quantities. Note first that 

dfl/dr = -H~ 3 dTi/dt + if I [1 - 12np + Antrf] (4.5) 

Next, a change of variables in (2.7) gives in the diagonal case: 

Tt=g u J fo(v 1 ,v 2 ,v 3 )(v l ) 2 (m 2 + g rs v r v s )- 1 / 2 (detg)- 1 / 2 dv 1 dv 2 dv 3 (4.6) 

Hence 

dTl/dt = (3k\ + kl + kfjTl 

+ 9 11 J fo(vi, V2, v 3 )(v 1 ) 2 F(v 1 , v 2 , v 3 )(m 2 + g rs v r v s )~ 3/2 (det g)~ 1/2 dv!dv 2 dv3 
where 



(4.7) 



F( Vl ,v 2 ,v 3 ) = (-g 11 ^) 2 - g 22 k 2 (v 2 ) 2 - g 33 kl(v 3 ) 2 ) (4.8) 

Note now that 

|^i^2,^3)| < (\k\\ + \k 2 \ + \k 3 \){m 2 + g rs v r v s ) (4.9) 

and so the integral in (4.7) can be bounded in modulus by 3HT^. Putting this information 
into (4.5) gives the desired bound. 

Lemma 4.4 Consider a maximally extended solution of (2.1)-(2.7) with m — 1 and / not 
identically zero. If gu — > 00 as r — > 00 then lim^oo T\j p = 0. If on the other hand on is 
bounded above and all gu are bounded below by a positive constant on an interval of the 
form [n, 00) then T\j p is bounded below by a positive constant on that interval. 
Proof It follows from (2.13) and (2.14) that T\ < Ca~ 2 p and this proves the first state- 
ment. To get the other conclusion, choose a cube C\ as in the proof of Lemma 4.2. Then 
T\ > Ca-^c' 1 while p < Ca-^c' 1 . Hence T\ j p > C> 0. 

5. The main results 

Lemma 5.1 (Compactness lemma) Let a sequence of reflection-symmetric global solutions 
of equations (2.1)-(2.7) be given. Then there exists a subsequence such that kj and p 
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converge uniformly on compact sets of R. T\ also converges uniformly on compact subsets 
(after possibly passing to a subsequence once more) and the limiting quantities satisfy 
(2.19) and (2.20). 

Proof The quantities kj are contained in the compact set K and so are, in particular, 
uniformly bounded. By (2.19) p is uniformly bounded. It follows that T* is uniformly 
bounded. Equation (2.20) now shows that d T kj is uniformly bounded. By Ascoli's theo- 
rem there exists a subsequence such that kj converges uniformly on the interval [—1,1]. 

Applying the theorem again shows that this subsequence has a subsequence such that M 
converges uniformly on [—2,2]. Continuing in this way we obtain a collection of subse- 
quences indexed by a positive integer n with the properties that for the nth subsequence 
kj converges uniformly on [— n, n] and each sequence is a subsequence of the previous one. 

The diagonal sequence has the property that kj converges on each compact subset of the 
real line. By the Hamiltonian constraint p also converges uniformly on compact sets along 
this subsequence. In the diagonal case the derivatives d r T\ are bounded, as was shown in 
Section 4 and applying Ascoli's theorem as before gives the remaining conclusions. 
Theorem 5.1 Let a global solution of equations (2.1)-(2.7) be given for which / is not 
identically zero. If at some time T\ the projection of the solution lies in the set Bi for some 
i then: 

(i) the projection lies in Bi for all r < t\ 

(ii) if there exists some T2 > T\ such that the projection of the solution lies in the comple- 
ment of Bi then it lies in the complement of all Bj for r > t^- 

(hi) as r — > — oo the projection converges to a point of the boundary which is not one of 
the points V~ 

Proof Suppose without loss of generality that i = 1. It follows from (2.20) that 

d T k\ = -Atx( P - trf )(3A; 1 1 - 1) - 87rtrTfci + 8^ (5.1) 

If at some time k\ < then the first and second terms on the right hand side of (5.1) 
are non-negative while the third term is positive. Hence d T k\ > 0. This implies the first 
conclusion of the theorem. Moreover it means that if the projection once leaves B\ it can 
never reenter it. A similar statement of course applies to any other Bj and this gives (ii). 
To prove (hi) note first that d T k\ is bounded below by a positive constant as long as p is. 
This shows that liminfT-^-oo p = 0. Equation (5.1) also implies that the integral of p on 
the interval (— oo, t±] must be finite so that for each % the integral of the right hand side 
of (2.20) is absolutely convergent. Hence each k\ tends to a limit as r — > — oo. By what 
has already been said it can only be a point of the boundary of K. The monotonicity of 
k\ shows that this limit cannot be one of the points V~ . 

Theorem 5.2 Let a global solution of equations (2.1)-(2.7) be given for which / is not 
identically zero. If at some time T\ the projection of the solution lies in the set A t for some 
i then as r decreases either: 

(i) the projection converges to a point of the boundary of K as r — > — oo or 

(ii) it reaches for some j or C at a finite time before T\ or 

(iii) it stays in Ai for all r < T\ and it has a point of one of the lines or the point C as 
an accumulation point. 
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Proof Suppose without loss of generality that i = 1. When the projection lies in Ai 
the inequalities k\ > k\, k\ > fcf and k\ > 1/3 hold. Suppose that on the time interval 
(—00, 7~i) the inequalities k\ — k% > A and k\ — k^ > A are satisfied for some A > 0. Then by 
Lemma 4.3 it follows that on this time interval T\jT^ and T\jT^ can be bounded below by 
a decreasing function which tends to 00 as r — > —00. Moreover, by Lemma 4.1 trT/p — > 1 
as r — > —00. Now define a sequence of solutions of (2.1)-(2.7) by it n (r) = u{r — n), 
r G (—00, n), where -u denotes any of the functions which make up the solution and n is 
a positive integer. By Lemma 5.1 there exists a subsequence such that kj and p converge 

uniformly on compact subsets. By the statements made above trT must tend to the same 
limit as p along this sequence. Also tends to this same limit and Tf and Tf tend to 
zero. Applying Lemma 5.1 again shows that the limits of these sequences satisfy (2.19) 
and (2.20). Because of the values of the limits they in fact satisfy the asymptotic system 
(3.1)-(3.3). The solution of the asymptotic system obtained inherits the properties that 
k\ > k\ and k\ > fcj. The only solutions of the aymptotic system which satisfy these 
inequalities on an interval of the form (— oo,7~i) are the vacuum solutions. If for some 
choice of subsequence this vacuum solution is not that corresponding to the point V{~ 
then the projection of the original solution must converge to that point, by Theorem 5.1. 
Otherwise every subsequence of the sequence of translated solutions has a subsequence 
which converges to the same solution of the asymptotic system. Hence the whole sequence 
converges to this solution and the projection of the original solution converges to V{~ . In 
both cases the solution of the original system converges a point of the boundary of K. 

It remains to consider the case where the above estimate is not satisfied for any A > 1. 
If the solution does not reach for some j or C in finite time then it stays in A t for all 
r < T\. Then it must have as an accumulation point either a point on for some j, C or 
V+ for some j. In the first two cases this gives case (iii) of the conclusion of the theorem. 
In the third case the solution enters Bj and so by Theorem 5.1 case (i) of the conclusion 
holds. 

In case (iii) of this theorem we can also consider a limit of translates of the solutions whose 
existence is guaranteed by Lemma 5.1. If the ratios a/b and a/c tended to zero as r — > — 00 
for the original solution then by Lemma 4.2 the ratios T^/Tf and T^/Tf would tend to 
infinity and the solution would belong to case (i). Thus in case (iii) it can be assumed 
without loss of generality (after possibly interchanging the indices 2 and 3) that b/a is 
bounded as r — > —00. Hence /^(^i —^2) i s finite. Since d T (k\ — fc|) is bounded it follows 

that k\ — k\ — > as r — > —00. Hence the solution obtained as a limit of translates has LRS 
geometry. It also satisfies p = trT. Information about the asymptotics of LRS solutions 
can thus be used to obtain information about the asymptotics of the solutions which fit 
into case (iii) of Theorem 5.2 but do not fit into case (i). 

Theorem 5.3 Let a solution of equations (2.19)-(2.20) be given which satifies the LRS 
condition k\ = k\. If at some time t\ the projection satisfies k\ < 1/3 then either: 

(i) the projection of the solution tends to the point ( — |, |, |) as r — > —00, 

(ii) it tends to the point (0, \) as r — > —00, 
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(iii) it reaches the point (|, |, g) at a finite time before T\ or 

(iv) it tends to the point (i, i, i) as r — > — oo 

Proof Suppose first that fcj < 1/3 — A for some A > 0. Then by Lemma 4.2 the ratio of 
T| = T| to increases without limit. Passing to a limit of translates in the familiar way- 
gives a solution of the equation d T k\ = —d>nk\p for which k\ satisfies the same inequality 
as before. There are only two such solutions, namely that for which p = and that for 
which k\ = 0. In the first of these cases the original solution must enter the region Bi and 
hence by Theorem 5.1 belong to case (i) of the conclusion of the present theorem. The 
only way of avoiding this is if, no matter which subsequence is chosen, the limiting value 
of k I is zero. Hence the projection of original solution must converge to the point (0, -|). 
Thus the solution belongs to case (ii) of the conclusion. Now consider the case where there 
is no A > with the given property. If, despite this, the ratio of T 2 2 to T± tends to infinity 
we can argue as before to show that the solution belongs to case (i) or (ii). If, on the other 
hand this ratio remains bounded then then ratio a/b must remain bounded and hence if 
k\ remains smaller than 1/3 for ever then J_^(l/3 — k\) is finite. It then follows as in 

the discussion following the proof of Theorem 5.2 that k\ — > 1/3. Thus the solution either 
belongs to case (iii) or case (iv). 

Theorem 5.4 Let a solution of the equations (2.1)-(2.7) with m = 1 and / not identically 
zero be given. Then k\ — > 1/3 and T\ / p — > for each i as r — > oo. 

Proof In [12] it was shown that the scale factors a, b and c are bounded below by a 
positive constant on any interval of the form [n,oo). Using (5.1) this statement can be 
strengthened. Suppose that k\ were negative on an interval of the form [ti,oo). Then 
it would follow, as in the proof of Theorem 5.1, that the integral of p on this interval 
was finite. But p is increasing on this interval, a contradiction. It follows that each k\ 
must become zero after a finite time and once this happens it must immediately become 
positive and stay positive. Thus, for n sufficiently large a, b and c are increasing. Consider 
now the behaviour of the quantity min{a, 6, c}. Suppose first that it tends to infinity as 
r — > oo. Then by Lemma 4.4 the ratios T\l p tend to zero as r — > oo. Construct a limit of 
translates as in the proof of Theorem 5.2 except that this time the translations should be 
done in the opposite direction. Then the limiting solution satisfies d T k\ = — Yl-Kp(k\ — |). 
This is the equation which is satisfied by a Bianchi I solution of the Einstein equations 
coupled to dust. It is well known and also easy to see directly that in the case of dust 
each k\ converges to 1/3 as r — > oo. Because a convergent subsequence can be extracted 
from any subsequence of the sequence of translates by integers it follows that k\ — > 1/3 
for the original solution of the Einstein-Vlasov system as well. Next consider the case 
where a is bounded on an interval of the form [ti,oo) while min{6, c} — > oo as r — > oo. 
Then by Lemma 4.4 the ratios Tf/p and Tf/p converge to zero as r — > oo while T\jp 
remains bounded away from zero. Equation (5.1) implies that d T k\ is bounded below by a 
positive constant if k\ < 1/3 and T^/T^ and Tg IT\ are less than 1 — A for some constant 
A > 0. However this contradicts the boundedness of a. Since the volume tends to infinity 
as r — > oo at least one of a, b or c must tend to infinity. It follows that to complete the 
proof we may assume without loss of generality that a and b are bounded while c tends to 
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infinity. By Lemma 4.4 T| / p — > while T^/p and Tf/p are bounded below by a positive 
constant. Now the integral fe|(r)dr is finite for z = 1, 2 and <9 T % is bounded. Hence 

k\ and fcf tend to zero as r — > oo and fcf — > 1. But the given behaviour of the pressures 
shows that for /cf > 1/3 and sufficiently large times d T k^ is negative, a contradiction. This 
completes the proof. 

6. A compact ificat ion 

This section is devoted to a finer examination of LRS solutions of the Einstein- Vlasov 
system with massless particles. For LRS solutions with k\ = k^ let k = k\, q = b/a, 
Q = Tl/ p. Then k\ = \(l - k) and p = (1/16tt)(± + k - f/c 2 ). Then the essential 
equations describing the dynamics are: 

d T k= l(l + 3k)(l-k)(Q-k) (6.1) 
d T q = l(l-3k)q (6.2) 

The quantity Q can be expressed entirely in terms of q and the initial data as follows: 

S fo(vi)vi(q 2 Vi + vi)~ 1/2 dv 1 dv 2 dv 3 ] 

Substituting (6.3) into (6.1) makes the equations (6.1) and (6.2) into an autonomous system 
of ordinary differential equations for k and q. Lemma 4.2 shows that Q(q) = 0(q 4 ^ 3 ) as 
q — > 0. This means that the system (6.1)-(6.2) can be extended in a C 1 manner to the 
boundary q = 0. Moreover, Q does not contribute to the linearization of the extended 
system at the critical point q = 0, k = 0. The eigenvectors of the linearization are directed 
along the k and q axes with eigenvalues —1/4 and 1/2 repectively. It follows (see e.g. 
[15]) that the dynamical system has an unstable manifold which is a curve tangent to the 
k axis at the point (0,0). This shows that for any initial value fo it is possible to find 
LRS solutions of the Einstein- Vlasov system where the distribution function has the initial 
value fo and where the quantities k\ converge to (0, |, |) as r — > — oo. Note that the stable 
manifold is just the k axis and so does not give rise to any smooth solutions of the Einstein- 
Vlasov system. The information about the linearization also determines the nature of the 
phase portrait near the singular point and shows that there are solutions for which k\ 
approaches zero but turns back before reaching it. A typical feature of Bianchi models is 
that the matter becomes dynamically negligeable near the singularity. No attempt will be 
made here to make this notion precise but one aspect of it is that the projection of the 
solution should tend to a point of the boundary of K as r — > — oo. The solution whose 
existence has just been shown is an exception to the rule. For anisotropic Bianchi I models 
with a perfect fluid no exceptional solutions of this kind exist. However they do occur for 
other Bianchi types[16]. 

The critical points of the system (6.1)-(6.2) in the region where < q and —1/3 < 
k < 1 are the points of the form (1/3, go) where qo has the property that Q(qo) = 1/3. 
Differentiating (6.3) and estimating the result in an elementary way leads to the inequality 
qQ' > Q(l — Q)- This shows that the function Q is strictly increasing for q > 0. Taking 



Q = q 2 
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account of the limiting values of Q, it follows that there is precisely one value qo for which 
Q(qo) = 1/3. Moreover, at this point qQ' > 2/9. The eigenvalues of the linearization at 

the corresponding critical point are — | ± \ \J\ — 2q Q'(q ). They both have negative real 
parts and the critical point is a sink. In particular no solution emerges from this critical 
point. Thus it is seen that if the projection of any LRS solution for massless particles 
approaches the point C as r — > — oo then the projection must stay at C for all time i.e. 
the solution must have isotropic geometry. This should be compared with the results of 
Newman [17] on isotropic singularities in solutions of the Einstein equations coupled to a 
radiation fluid. 

Appendix: Fluids with nonlinear equation of state 

Consider a perfect fluid with equation of state p = f(p) which satisfies the following 
general assumptions: 

(i) / is a continuous function from [0, oo) to itself with /(0) = which is C 1 for p > 

(ii) < f'(p) < 1 for all p > 

(iii) there exists a constant C < 1 such that p < Cp for p < 1 

Assumptions (i) and (ii) are standard. The third assumption is, when (i) and (ii) are 
satisfied, equivalent to the assumption made in [7] that the solution is not asymptotically 
stiff at low densities. In the case of a linear equation of state f(p) = kp the assumptions 
(i)-(iii) are satisfied if and only if < k < 1. In a Bianchi I spacetime it follows from 
the momentum constraint (2.2) that the four-velocity of the fluid is orthogonal to the 
hypersurfaces of homogeneity. Hence the energy density p measured by an observer whose 
wordline is orthogonal to the hypersurfaces of homogeneity is the same as that measured by 
a comoving observer. Equations (2.19) and (2.20) are valid as in the case of the Einstein- 
Vlasov system. For a fluid it can be assumed without loss of generality that the solution is 
reflection-symmetric because given any initial data, it suffices to do a linear transformation 
of the coordinates which simultaneously diagonalizes the metric and second fundamental 
form in order to transform the given data to data for a reflection-symmetric spacetime. 

In the case of a fluid T- = pS* and hence T- = pdp where p = p/H 2 . If the equation of 
state is linear then p can be expressed as a function of p alone and (2.19) and (2.20) then 
reduce to a system of ordinary differential equations which suffice to determine kj from 
initial data. For a nonlinear equation of state this is no longer the case. The equations 
(2.19) and (2.20) no longer form a closed system and must instead be considered as the 
projection of a bigger system, as in the case of the Vlasov equation. This is one reason 
why the linear case has been studied preferentially in the literature. Nevertheless it turns 
out that the projection can be analysed very effectively in the general case. 

The first question which needs to be addressed is that of global existence in r i.e. the 
equivalent of Lemma 2.1 for a fluid. This follows from the results of [7]. The assumption 
(iii) has been used at this stage. A direct calculation shows that the quantity (k\ — 
k?,) I (k\ — /eg) is independent of r whenever k\ — /cf is non-zero. Moreover if k\ — /cf is zero 
at some time it remains zero. Hence the projection of each solution is constrained to move 
on a straight line in K passing through the centre C. This is already a much stronger 
statement than could be proved in the case of the Vlasov equation. To find out how the 



15 



projection moves on this straight line, the time derivative of the dimensionless version of 
the density will be calculated. For a fluid it is given by 

d T p = (p-p)(l-24np) (Al) 

Noting that the Hamiltonian constraint implies that 2Aixp < 1 with equality only at the 
point C, it can now be seen that the projection moves monotonically from the boundary of 
K at r = — oo to the centre C of K at r = oo. This qualitative behaviour is independent of 
the equation of state satisfying (i)-(iii). The only difference is in the speed with which the 
projection moves along the radial line at different times. Equation (Al) also makes clear 
that this picture changes completely if the equations of state considered here are replaced 
by the limiting case of a stiff fluid, p = p. 
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